Algorithm of Pulmonary Vascular Segment and Centerline Extraction

Because pulmonary vascular lesions are harmful to the human body and difficult to detect, computer-assisted diagnosis of pulmonary blood vessels has become the focus and difficulty of the current research. An algorithm of pulmonary vascular segment and centerline extraction which is consistent with the physician's diagnosis process is proposed for the first time. We construct the projection of maximum density, restore the vascular space information, and correct random walk algorithm to satisfy automatic and accurate segmentation of blood vessels. Construct a local 3D model to restrain Hessian matrix when extracting centerline. In order to assist the physician to make a correct diagnosis and verify the effectiveness of the algorithm, we proposed a visual expansion model. According to the 420 high-resolution CT data of lung blood vessels labeled by physicians, the accuracy of segmentation algorithm AOM reached 93%, and the processing speed was 0.05 s/frame, which achieved the clinical application standards.


Introduction
In recent years, the air is in poor quality, which seriously affects people's life and health, resulting in an increasing incidence of lung cancer. Because the early symptoms of lung cancer are relatively minor, it is not easy to be discovered, and 80% of lung cancer patients are already in the middle and late clinical stage when they are confirmed, missing the best period of radical operation. If the lung cancer can be detected earlier and get differential diagnosis by category, it can receive standard treatment before the lesion spread, the patient's five-year survival rate can reach more than 60% [1,2]. If the lung CT data which intuitively reflect the lung condition can be accurately analyzed, lung disease can also be prevented and diagnosed, and more human life can be saved. Pulmonary vascular lesions are one of the early signs of lung cancer, so accurate division and observation of pul-monary vascular areas are very important for diagnosis and treatment [3].
Due to the small morphological structure of pulmonary blood vessels and complex distribution environment, the process of imaging is easily affected by noise and surrounding tissues. This has brought great difficulty to the complete extraction of lung blood vessels, which is also the focus of the current research. According to the characteristics of the vascular section obeying Gaussian distribution, Foncubierta-Rodríguez et al. [4] constructed an enhancement model of vascular segmentation to determine the direction of vascular growth and extract blood vessels. According to the vascular anatomy, Orkisz et al. [5] calculated the global threshold to divide the lung image, and then, the blood vessels are extracted by discriminating the connectivity. Based on the two-dimensional CT image, Ibrahim et al. [6] selected the initial point manually, calculated the blood vessel position step by step, and then segmented the blood vessel. Lai et al. [7] enhanced the vascular response by filtering the image from a 3D angle. Gong et al. [8] extracted the blood vessels by constructing the selfadaptive threshold according to the difference between the gray level of the blood vessel and the surrounding area. Charbonnier et al. [9] constructed a matching model in accordance with blood vessels extending throughout the lungs to segment the blood vessels. Based on the local feature model, Gupta et al. [10] predicted the blood vessel direction. According to the samples of axial position, coronal position, and sagittal position, Phellan et al. [11] constructed a CNN model to segment the vascular area. Lidayová et al. [12] established the model of normal blood vessel and abnormal blood vessel to realize the segmentation of blood vessel. Jawaid et al. [13] improved level set to extract the vascular boundary information. Xiao et al. [14] used the grayscale and shape characteristics to constrain and segment blood vessels. Flórez and Orkisz [15] established the cylindrical model to simulate the distribution of blood vessels. Hu et al. [16] improved the crossaccess law function to extract blood vessels.
Because it is difficult for doctors to observe small changes of blood vessels through the naked eye to make accurate diagnosis, computer-assisted diagnosis is required. Based on the above extraction of the vascular centerline, the mainstream diagnostic method is as follows: adjust the window width and window position [17] when blood vessels are observed, the 3D model [18] shows the overall external morphology of blood vessels, a virtual endoscope [19] shows the inside of the blood vessels. The above algorithms cannot give the internal and external characteristics of blood vessels at the same time, which limited the ability of helping physician to diagnose.
The main problems of vascular division, centerline extraction, and auxiliary diagnosis can be summarized as follows: (1) the 3D model of extracting vascular segment is complex, and it needs a large amount calculation, so it cannot be calculated and displayed in real time; (2) it is difficult to divide the vascular area accurately just by gray scale; and (3) there is no accurate and intuitive way to show the internal and external conditions of blood vessels.
For those reasons, we study the process of doctors' diagnosis of vascular lesions and construct a diagnostic model for doctors' vision diagnosis according to the theory of medicine, anatomy, and image graphics: (1) establish a maximum density projection model based on vascular anatomy and imaging principles and significantly reduce the amount of data while retaining the area where the blood vessels are located; (2) adaptive random walk algorithm is constructed by combining gray level with local information; and (3) construct the mechanism of vascular centerline extraction and vascular expansion and display and diagnose vascular from multiangle.

Details of Algorithm
When physicians detect pulmonary vascular lesions, removing the interference of the examination bed, muscles, and soft tissues is the first task. They focus on the lung parenchyma and then determine the starting and ending points of the blood vessels to be observed based on the local images of limited frames, achieving the extraction of blood vessels. Finally, they extract the vascular centerline and give diagnosis based on the anatomical knowledge.
The computer simulates how doctors diagnose lesions as follows: the lung image sequence was firstly pretreated to remove noise and check a bed, and the image was simplified by extracting the lung parenchyma area. Then, it determined the starting and ending points of blood vessels to be observed using global and local models to constrain them and took the characteristics of anatomy and image graphics into consideration to realize the extraction of blood vessels' central line. Figure 1 shows the flow of pulmonary vascular CAD segmentation.
2.1. Pretreatment. CT images can intuitively show the morphological structure of human tissues (blood vessels, musculoskeletal, etc.). According to the physician's diagnosis process, the blood vessels that physician concerned only exist in the lung parenchyma. The ends of the blood vessels are tiny; they are only displayed as a limited number of pixels in the image and are vulnerable to noise. Therefore, the image needs to be preprocessed.
2.1.1. Noise Reduction. CT image always includes pepper and salt noise [20] which presented as a bright discrete pixel; it interferes with the extraction of pulmonary blood vessels' tip and needs to be removed. According to the characteristics that noises are distributed as discrete points, a twodimensional median filtering algorithm can be used to remove them. The formula is as follows: where f ðx, yÞ and gðx, yÞ are the gray values of the pixels with coordinates (x, y) in the original image and the processed image, respectively; W is a two-dimensional filter template; and med is a median filter function.

Lung Parenchyma Extraction.
The area of blood vessels that doctors care about is inside the lung parenchyma. Whether the extraction of the blood vessels is good or bad directly affects the subsequent diagnosis. Therefore, the computer needs to focus on the area where the pulmonary nodules are located. By analyzing the statistical distribution of lung image pixels, Qiu et al. [21] constructed an optimized Gaussian dual-mixing model to calculate the global segmentation threshold quickly, and the complete lung parenchyma region can be extracted. In order to meet the clinical application standards of accuracy and processing speed, the algorithm in literature [21] is applied to extract the lung parenchyma region and confirm the left and right lung regions.

Set the Starting and Ending Points of Blood Vessels.
When the doctor identifies the focal point, he will determine the observed vascular area through clinical needs. When the computer simulates this process, according to the morphological structure of blood vessels, it sets the blood vessel's starting point P b and the end point P e whose corresponding coordinates are (x b , y b , z b ) and (x e , y e , z e ) on the computer image through human-computer interaction, and z b ≤ z e .

Rough Extraction by a Local 3D Model.
The physician only cares about the blood vessel between point P b and point P e , so the blood vessel area needs to be extracted.
Since blood vessels are of connectivity, as shown in Figure 2, only handle ½z b , z b + 1, ⋯z e is not rigorous, where m is the point on the blood vessel.
Anatomically, blood vessels are spatially continuous. From the perspective of image imaging, blood vessels appear as bright and isolated circular areas in two-dimensional CT images.
In order to restore the anatomical characteristics of blood vessels, we introduce the maximum density projection which is integrated with the local 3D information to greatly reduce the amount of data while ensuring the information of blood vessels.
The maximum intensity projection is a method for 3D data that projects in the visualization plane the voxels with maximum intensity that fall in the way of parallel rays traced from the viewpoint to the plane of projection. The equation is as follows: where MIPnðx, yÞ is the grayscale value at the midpoint (x, y ) of the nth frame of the MIP image; H and W are the horizontal and vertical resolutions, respectively; SN r is the number of projected layers; and Nm is the sum of original CT layers. I k ðx, yÞ is the grayscale value at the point (x, y) in the kth layer of the original CT sequence images. Figure 3 is a schematic diagram of the MIP projection of the twodimensional CT data. The blood vessels in the MIP image are locally highlighted and continuous areas, as shown in Figure 4(a). Subsequent extraction of vascular areas is needed to constrain the vascular centerline method and reduce the amount of data. In the MIP image, P b and P e are chosen as seed points, The area of the blood vessel is obtained by using C -means clustering [22,23], and then, the image of the suspected area I i as Figure 4(b) is obtained through the 8 connected areas. Follow-up study is based on I i .

Segment by Automatic Random Walk
Algorithm. Random walk (RW) algorithm [24] regarded the twodimensional image as a connected undirected weighted graph containing fixed vertices and edges. The unlabeled pixels begin to walk along the edge from the vertices. According to the maximum probability of each pixel arriving at each labeled pixel, the class of each vertex is judged.
The main steps are (1) manually mark the seed point; (2) establish a random walk model G = ðV, EÞ, where V is the set of vertices of the image and E is the set of undirected   In recent years, Dong et.al [25] introduced subMarkov to realize image segmentation. Zheng et al. [26] constrain RW according to target characteristics. Bui et al. [27] established a 3D model to realize the heart division. All of the above methods get good results. As blood vessels run through multiple layers of images, it is inefficient to manually mark seed points on each layer of images. In a 2D image, blood vessels are shown as limited sections of pixels, and the number of pixels is small, so the accuracy of seed selection will directly affect the segmentation effect. For this reason, we study the selection of seed and weight calculation and then propose an algorithm.
2.3.1. Adaptive 3D Seed Selection. The traditional RW algorithm segments the image semiautomatically. The user needs to manually mark the target seed and background seed on the 2D image, and then, the probability of each vertex to these two kinds of seed is used to judge which one they belong to, thus dividing the image into the target area and the background area.
The blood vessel exists on the CT image as a circular or cylindrical section. When using the traditional RW algorithm, the physician needs to spend a lot of time on accurately marking the seed on each 2D image. Besides, the number of blood vessel pixels is limited; the physician may mark the wrong seed because of a slight deviation, which cannot guarantee the completeness of the segmentation.
According to the characteristics of angiology, anatomy, and CT imaging, we optimize the process of selecting seeds as Figure 5. From the medical point of view, the components contained in the blood vessels are uniform, so the pixel values are similar in CT image. Anatomically, it is a continuous region in a 3D space. From the aspect of imaging, blood vessels are locally bright area.
The approximate area of the blood vessel in the 3D space is marked as I i , with the point P i (P 0 = P b ) as the target seed and the nonvascular region in I i as the background seed, and the divided blood vessels of the current layer was recorded as C i . According to the continuity of vascular, C i ∩ C i+1 ≠ ∅, according to the uniformity of blood composition, the pixel value of the blood vessel should change a little; then, the coordinate of the central point of C i ∩ C i+1 is assigned as P i+1 .

Calculate Weights.
The weights determine that the probability of each vertex belongs to different classes. The traditional RW algorithm only considers the grayscale information between pixels and does not consider the geometric information between pixels. For more accurate segmentation, the improved weight function is    Computational and Mathematical Methods in Medicine where gð:Þ represents the grayscale, hð:Þ represents the coordinate value, α is the weight of grayscale, ð1 − αÞ is the weight of distance, and L i,j is the matrix of distance. RW algorithm is rebuilt to realize the automatic segmentation of blood vessels. The selected seeds are few, the calculation area is small, and the speed is fast. It reduces the workload of doctors and improves the diagnostic efficiency.

Centerline Extraction and Expansion.
After the doctor confirms the vascular area, he reconstructed the 3D characteristics of the blood vessels and set up the center point to observe the pathological conditions of the blood vessels. Because the spacing of CT data layers and the image resolution are not uniform, in order to restore the true situation of the lungs, the computer needs to interpolate the data to develop the isotropic data; that is, the lung data observed from any direction are uniform. According to MIP image processing, it can be seen that the space of blood vessels is limited, so it needs to be interpolated. The gray scale of the point P i is used as the seed to cluster on the corresponding layer of the original image sequence to obtain the blood vessel area, whose center and the inscribed circle's radius R are calculated. Then, the point P i in the image is seen as the center of the ball with a radius of 2R; the data sequence is interpolated vertically according to the resolution of the axis to form isotropic data.
Calculate the three eigenvalues λ 1 , λ 2 , and λ 3 and the eigenvectors ν 1 , ν 2 , and ν 3 , which represent the trend of movement in different directions, as shown in Figure 6. Take ν 1 as the blood vessel direction and 2R as the search radius; make some adjustment R to find the position of next blood vessel to obtain a complete blood vessel area.
In the process of the calculation, the following four situations occur, as shown in Figure 7, in which the red strip is a blood vessel and the balls show the adjustment process of centerline.
Case 1. Start off from P b ; calculate ν 1 . In the direction of ν 1 , jump with the step 2R 0 . Calculate the coordinates of the center of ball's gravity P 1 and radius R 1 and save it.
Case 2. Start off from P 1 ; calculate ν 1 . In the direction of ν 1 , jump to P 1 ′ with the step 2R 1 . Constantly adjust the center of gravity and radius until they do not change; save the center of gravity P 2 and radius R 2 .
Case 3. Start off from P 2 ; calculate ν 1 . In the direction of ν 1 , jump to P 2 ′ with the step 2R 2 , and the circle centered on P 2 ′ with the radius R 2 which is not in the MIP image, indicating that the step length is too large, which needs to be adjusted. When the cross section of circles reaches 10%, the step length is set to R 2 ′. The coordinates at this point are P 2 ″. Adjust the center of gravity and radius continuously until they do not change; save the center of gravity P 3 and radius R 3 .
Case 4. Start off from P 3 ; calculate ν 1 . In the direction of ν 1 , jump to P 3 ′ with the step 2R 3 . If P e is exceeded, take P e as the center of the circle and adjust the length of the step, and the center of gravity of the ball is calculated and saved.
When the physician determines the vascular centerline, he needs to observe the internal and external conditions of the blood vessel; this process is simulated by the computer as follows: first, the blood vessels are straightened according to the centerline, and then, according to the anatomical principles, the blood vessels are expanded according to 0°4 5°13 5°P b P 1 P 2 P 3 P 4 P 5 P 6 P 7 P 8 P e 90°P b P 1 P 2 P 3 P 4 P 5 P 6 P 7 P 8 P e

Experiments and Result Analysis
The experiment used 50 sets of CT data with different resolutions and thickness collected in 2020; 420 blood vessel areas were extracted. The programs were compiled in VS2015 with the WIN7 operating system.
3.1. Verify the Effect of MIP Algorithm. The same group of CT sequences was reconstructed as the layer thickness is 3 mm, 9 mm, and the entire sequence, respectively, by the MIP algorithm and compared with the original sequence. The effect is shown in Figure 9. The MIP algorithm can restore the strip shape of vascular and reduce the risk of wrong detection. With the increase of layer thickness, it can display the distribution of peripheral blood vessels more clearly and ensure that the follow-up extraction of vascular centerlines is more complete.   To verify the effect of blood vessel extraction, three experts who have been engaged in medical imaging for many years were asked to mark blood vessels to verify the effect of the algorithm from the quantitative and qualitative levels.
The quantitative aspects are as follows: take area overlap measure (AOM) [28] as an evaluation indicator of the seg-mentation effect; it is defined as where AOM is the area overlap measure, A is the image marked by the doctor, B is the image segmented by   Table 1.
The qualitative aspects are as follows: using medical diagnostic marks: "accurate" means that the vascular division is complete and almost identical with the area marked by experts; "general" indicates that there is a some deviation between the results of the two, but it does not affect the diagnosis; and "poor" indicates that there is a large deviation between the results of the two, which affects the diagnosis. The results are shown in Table 2.
This paper uses the Contrast-Enhanced (CE) algorithm [4] from both qualitative and quantitative aspects to construct a vascular enhancement model and achieve vascular segmentation. 3D algorithm [7] established a 3D model from the inside of the blood vessels. Multidirectional (MD) algorithm [11] observed the blood vessels from the axial, coronal, and sagittal direction, determined the direction of blood vessels, and then segmented them. Level set (LS) algorithm [13] transformed the segmentation problem into the problem of internal and external forces' balance. RW [24] selected the seed points manually on each layer of images to achieve vascular segmentation, not considering the time spent. Our algorithm takes full account of vascular anatomy and imaging principles, proposes an automatic selection algorithm for seed points, and optimizes the weight function; although a little more time was spent than CE [4] algorithm, the segmentation effect AOM reached 93%, and the effect that achieves the "general" was the best, thus confirming the effectiveness of the proposed algorithm.
To demonstrate the effect of our algorithm, we extract the whole pulmonary blood vessels. It is very difficult. The algorithm proposed in this paper can extract most of the blood vessels, but there are still some incomplete blood vessels extracted. As shown in Figure 10(b), the blood vessels are not completely extracted in the area where they enter the lungs. The main reason is that the left and right lungs of this sequence are relatively close. During the CT imaging process, the left and right lungs are not effectively divided. Further research is needed on this issue.

The Extraction of Vascular Centerline.
In order to verify the effect of the extraction of vascular centerline and show the observation effect of the current major assisted diagnosis method, our algorithm, 3D display, and virtual insight are compared. The 3D display is shown in Figure 11(a); only the external condition of the blood vessel can be observed; the virtual insight is shown in Figure 11(b); only the internal condition of the blood vessel can be observed; and our algorithm is shown in Figure 11(c), which shows the inside and outside of the blood vessel at the same time and intuitively shows the effect of extraction of vascular centerline.
In order to intuitively demonstrate the effect of vascular expansion, we chose unbranched blood vessels, branching blood vessels, and blocked blood vessels, which were displayed at 0°, 45°, 90°, and 135°. Unbranched blood vessels, as shown in Figure 12, are spread in different directions; the centerline is in the center of the blood vessel. For branch-ing blood vessels, as shown in Figure 13, the longest blood vessel (red) and the shortest blood vessel (blue) are selected for straightening. It can be seen that this algorithm is not affected by the blood vessels next to the branch when extracting the central line. The central line can be extracted accurately. Figures 12 and 13 show the effects of normal vascular centerline's extraction and spread, and the interior of blood vessel is relatively smooth. Figure 14 shows the  Computational and Mathematical Methods in Medicine spreading effect of blocked blood vessel. It can be seen that the algorithm can accurately extract the centerline, and a lot of blocked spots can be observed inside the blood vessels. From the above demonstration, the algorithm presented in this paper can accurately extract the vascular centerline and directly display the internal and external conditions of blood vessels, so it can assist the physician to make an accurate diagnosis.

Conclusion
Computer-aided detection of pulmonary blood vessels has become a hot topic and difficulty in research. In this paper, the segmentation algorithm of pulmonary vascular CT image is proposed, which combines local 3D information and enhances vascular area. Optimize the random walk algorithm to meet the requirements of separating blood vessel accurately. We applied the knowledge of pathology, anatomy, and image graphics to the vascular centerline extraction and proposed the method of spreading blood vessels which can show the internal and external characteristics of blood vessels more clearly, thus assisting physicians in making accurate diagnoses. It is of great value on clinical application and lays a good foundation for the follow-up detection of vascular lesions. However, the vascular division is incomplete when the left and right lungs are near, and further research is needed.

Conflicts of Interest
The authors declare that they have no conflicts of interest.